!-------------------------------------------------------------LICENSE--------------------------------------------------------------!
!                                                                                                                                  !
!The MAP code is written in Fortran language for magnetohydrodynamics (MHD) calculation with the adaptive mesh refinement (AMR)    !
!and Message Passing Interface (MPI) parallelization.                                                                              !
!                                                                                                                                  !
!Copyright (C) 2012                                                                                                                !
!Ronglin Jiang                                                                                                                     !
!rljiang@ssc.net.cn                                                                                                                !
!585 Guoshoujing Road. Pudong, Shanghai, P.R.C. 201203                                                                             !
!                                                                                                                                  !
!This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License         !
!as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.             !
!                                                                                                                                  !
!This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of    !
!MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more details.                        !
!                                                                                                                                  !
!You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software     !
!Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.                                                   !
!                                                                                                                                  !
!-------------------------------------------------------------LICENSE--------------------------------------------------------------!

!==================================================================================================================================|
subroutine extremum_calculate(var_min, var_max, ro, gp, mp, te, beta, vx, vy, vz, bx, by, bz, jx, jy, jz,                          &
               fgx, fgy, fgz, fggx, fggy, fggz, fmgx, fmgy, fmgz, fmtx, fmty, fmtz, var_flag, n_var, nx, ny, nz, ng, n_var_total)
!==================================================================================================================================|

   implicit none

   include 'data_manage_precision.include'

   integer(4), intent(in) :: n_var, nx, ny, nz, ng, n_var_total
   logical(4), intent(in) :: var_flag(n_var_total)
   real(output_real8_flag), intent(inout) :: var_min(n_var), var_max(n_var)
   real(output_real8_flag), dimension(nx, ny, nz), intent(in) :: ro, gp, mp, te, beta, vx, vy, vz, bx, by, bz,                     &
      jx, jy, jz, fgx, fgy, fgz, fggx, fggy, fggz, fmgx, fmgy, fmgz, fmtx, fmty, fmtz

   integer(4) :: i, j, k, i_var, j_var

!----------------------------------------------------------------------------------------------------------------------------------|
!  ro, gp, mp, te, beta, vx, vy, vz, vv, bx, by, bz, bv, jx, jy, jz, jv, fgx, fgy, fgz, fgv, fggx, fggy, fggz, fggv, 
!  fmgx, fmgy, fmgz, fmgv, fmtx, fmty, fmtz, fmtv
!----------------------------------------------------------------------------------------------------------------------------------|
   do k = ng, nz - ng
   do j = ng, ny - ng
   do i = ng, nx - ng
      i_var = 1
      j_var = 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), ro(i, j, k))
         var_max(i_var) = max(var_max(i_var), ro(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), gp(i, j, k))
         var_max(i_var) = max(var_max(i_var), gp(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), mp(i, j, k))
         var_max(i_var) = max(var_max(i_var), mp(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), te(i, j, k))
         var_max(i_var) = max(var_max(i_var), te(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), beta(i, j, k))
         var_max(i_var) = max(var_max(i_var), beta(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), vx(i, j, k))
         var_max(i_var) = max(var_max(i_var), vx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), vy(i, j, k))
         var_max(i_var) = max(var_max(i_var), vy(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), vz(i, j, k))
         var_max(i_var) = max(var_max(i_var), vz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(vx(i, j, k) * vx(i, j, k) +                                                     &
                                           vy(i, j, k) * vy(i, j, k) +                                                             &
                                           vz(i, j, k) * vz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(vx(i, j, k) * vx(i, j, k) +                                                     &
                                           vy(i, j, k) * vy(i, j, k) +                                                             &
                                           vz(i, j, k) * vz(i, j, k)))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), bx(i, j, k))
         var_max(i_var) = max(var_max(i_var), bx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), by(i, j, k))
         var_max(i_var) = max(var_max(i_var), by(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), bz(i, j, k))
         var_max(i_var) = max(var_max(i_var), bz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(bx(i, j, k) * bx(i, j, k) +                                                     &
                                             by(i, j, k) * by(i, j, k) +                                                           &
                                             bz(i, j, k) * bz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(bx(i, j, k) * bx(i, j, k) +                                                     &
                                             by(i, j, k) * by(i, j, k) +                                                           &
                                             bz(i, j, k) * bz(i, j, k)))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), jx(i, j, k))
         var_max(i_var) = max(var_max(i_var), jx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), jy(i, j, k))
         var_max(i_var) = max(var_max(i_var), jy(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), jz(i, j, k))
         var_max(i_var) = max(var_max(i_var), jz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(jx(i, j, k) * jx(i, j, k) +                                                     &
                                             jy(i, j, k) * jy(i, j, k) +                                                           &
                                             jz(i, j, k) * jz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(jx(i, j, k) * jx(i, j, k) +                                                     &
                                             jy(i, j, k) * jy(i, j, k) +                                                           &
                                             jz(i, j, k) * jz(i, j, k)))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fgx(i, j, k))
         var_max(i_var) = max(var_max(i_var), fgx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fgy(i, j, k))
         var_max(i_var) = max(var_max(i_var), fgy(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fgz(i, j, k))
         var_max(i_var) = max(var_max(i_var), fgz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(fgx(i, j, k) * fgx(i, j, k) +                                                   &
                                             fgy(i, j, k) * fgy(i, j, k) +                                                         &
                                             fgz(i, j, k) * fgz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(fgx(i, j, k) * fgx(i, j, k) +                                                   &
                                             fgy(i, j, k) * fgy(i, j, k) +                                                         &
                                             fgz(i, j, k) * fgz(i, j, k)))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fggx(i, j, k))
         var_max(i_var) = max(var_max(i_var), fggx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fggy(i, j, k))
         var_max(i_var) = max(var_max(i_var), fggy(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fggz(i, j, k))
         var_max(i_var) = max(var_max(i_var), fggz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(fggx(i, j, k) * fggx(i, j, k) +                                                 &
                                             fggy(i, j, k) * fggy(i, j, k) +                                                       &
                                             fggz(i, j, k) * fggz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(fggx(i, j, k) * fggx(i, j, k) +                                                 &
                                             fggy(i, j, k) * fggy(i, j, k) +                                                       &
                                             fggz(i, j, k) * fggz(i, j, k)))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fmgx(i, j, k))
         var_max(i_var) = max(var_max(i_var), fmgx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fmgy(i, j, k))
         var_max(i_var) = max(var_max(i_var), fmgy(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fmgz(i, j, k))
         var_max(i_var) = max(var_max(i_var), fmgz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(fmgx(i, j, k) * fmgx(i, j, k) +                                                 &
                                             fmgy(i, j, k) * fmgy(i, j, k) +                                                       &
                                             fmgz(i, j, k) * fmgz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(fmgx(i, j, k) * fmgx(i, j, k) +                                                 &
                                             fmgy(i, j, k) * fmgy(i, j, k) +                                                       &
                                             fmgz(i, j, k) * fmgz(i, j, k)))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fmtx(i, j, k))
         var_max(i_var) = max(var_max(i_var), fmtx(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fmty(i, j, k))
         var_max(i_var) = max(var_max(i_var), fmty(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), fmtz(i, j, k))
         var_max(i_var) = max(var_max(i_var), fmtz(i, j, k))
         i_var = i_var + 1
      endif
      j_var = j_var + 1
      if (var_flag(j_var)) then
         var_min(i_var) = min(var_min(i_var), sqrt(fmtx(i, j, k) * fmtx(i, j, k) +                                                 &
                                             fmty(i, j, k) * fmty(i, j, k) +                                                       &
                                             fmtz(i, j, k) * fmtz(i, j, k)))
         var_max(i_var) = max(var_max(i_var), sqrt(fmtx(i, j, k) * fmtx(i, j, k) +                                                 &
                                             fmty(i, j, k) * fmty(i, j, k) +                                                       &
                                             fmtz(i, j, k) * fmtz(i, j, k)))
      endif
   enddo
   enddo
   enddo

!----------------------------------------------------------------------------------------------------------------------------------|
   return
end subroutine extremum_calculate
